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Abstract 

Nonlinear time-dependent differential equations for the Hele-Shaw, Saffman- 
Taylor problem are derived. The equations are obtained using a separable 
ansatz expansion for the stream function of the displaced fluid obeying a Dar- 
cian flow. Suitable boundary conditions on the stream function, provide a 
potential term for the nonlinear equation. The limits for the finger widths de- 
rived from the potential and boundary conditions are 1 > A > in units of 
half the width of the Hele-Shaw cell, in accordance with observation. Station- 
ary solutions with no free phenomenological parameters are found numerically. 
The dependence of asymptotic finger width on the physical parameters of the 
cell compares satisfactorily with experiment. The correct dispersion relation 
for the instabilities is obtained from the time dependent equation. 

PACS numbers: 47.20.Dr, 47.54.+r, 68.10.-m 
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1 Introduction 

The Hele-Shaw cell experiment of a less viscous fluid displacing by a more viscous one 
is the paradigm of fingering phenomena [|l|. The time-honored, work of Saffman and 
Taylor 10] showed both experimentally and theoretically that fingers with beautifully 
regular profiles arise in the Hele-Shaw cell. 

In the intervening years a large body of works has added to our knowledge of the 
fingering phenomenon in various branches of the dynamics of continuous media, such 
as dendritic growth, directional solidification, diffusion- initiated aggregation, flame 
front propagation, electromigration, as well as fingering in porous media. The porous 
medium problem prompted the initial research of fingering phenomena. The topic is 
more than ever now, of the utmost importance in problems of transport in saturated 
and unsaturated soils, groundwater pollution, etc. [] Moreover, research in the field is 
now evolving into more complicated geometries, like that of a spherical shape cell. Q 

Despite considerable efforts, the fingering phenomenon remains in many aspects 
uncharted territory. As Tanveer^ described it very lucidly, the phenomenon defies 
intuition. In the words of Tanveer]^ , the disproportionally large influence of small 
effects, like local inhomogeneities, thin film effects, and the crucial role of surface 
tension make the theoretical description very difficult. So crucial is the surface ten- 
sion in the fingering phenomenon, that the problem is ill-posed without it0. The 
original work of Saffman and Taylor p[ observed fingers in a cell with a less viscous 



^The website http: / /www.maths.ox.ac.uk^howison/Hele-Shaw /helearticles.bib , cited in ref. 



carries an extensive (more than 600 papers), up to date list of references on the Hele-Shaw problem. 



fluid (oil, air, water), penetrating into a more viscous fluid (oil, glycerine). An initial 
instability develops into a finger or several competing fingers. Saffman and Taylor 0] 
compared the experimentally photographed finger profiles to a theoretical prediction. 
The formula was obtained by taking advantage of the analytical properties of the 
complex fluid potential (potential and stream function). Reasonable boundary condi- 
tions, and neglect of the surface tension, gave good agreement with the experimental 
profiles only in the limit of very large capillary number (or alternatively small sur- 
face tension). A drawback of this theoretical approach lies in the absence of surface 
tension that, causes mathematical ill-posedness of the problem. This fact is reflected 
in the absence of a limitation for the asymptotic size of the finger. Contrarily to the 
measured profiles, that were always found to be limited by approximately one half of 
the cell width, the theoretical expression showed no such lower bound. 

McLean and Saffman]^ addressed the problem including surface tension. Simi- 
larly to Saffman and Taylor they emphasized the role of analyticity. Tanveer0 has 
shown that the approach of McLean and Saffman|]^, is equivalent to an expansion 
in a parameter related to the finger width, capillary number and aspect ratio. The 
expansion fails at the tail of the finger and higher order terms are needed. The re- 
sults of McLean and Saffman |]^ nevertheless yielded profiles that matched very well 
the front (nose) part of the finger. Also they found a limit of A > | to the finger 
width. Their curve for the dependence of the finger width on the capillary number 
show the right trend, albeit below the measured values. There appear to be, however, 
other branches of solutions (two at least) @], that compete with the branch found by 
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McLean and Saffman|0]. These branches are presumably unstable as compared to 
the one found by McLean and Saffman0. 

More accurate experiments with many different aspect ratios ( the ratio of width 
of the channel with respect to its thickness) showed that the limit is more likely 
around A > 0.45 of the cell widthp]. 

Amongst the phenomenological approaches to the problem, the remarkable work 



of Pitts|TO|, took advantage of the observed dependence of the curvature of the finger 
on the angle, to obtain an analytical scale covariant expression for the finger shape 
that fitted measured values (by Saffman and Taylor and by Pitts himself) extremely 
well, especially for finger widths smaller than A ~ 0.8. For wider fingers, the profile 
function was found to miss the measured finger by a small amount. To this day, there 
is no other stationary finger analytical expression in the literature that fits the data. 
Pitts derived also an expression for the dependence of the finger's asymptotic width 
on the capillary number that includes a phenomenological (fitted) parameter taking 
into account the fluid films left behind by the passing finger. The fit was also very 
successful. The theoretical basis for the phenomenological assumptions used by Pitts 
are nevertheless unclear. 



DeGregoria and Schwartz ||11|] used a boundary integral Cauchy technique to inves- 



tigate the production and propagation of fingers in the Hele-Shaw cell. Tryggvason 



and Aref ||T2[ used a vortex in line method to determine the relationship between fin- 
ger width and flow parameters. Both works, as well as the stationary calculation of 
McLean and Saffman give similar results concerning the finger width dependence 



on capillary number. A marked improvement is found in the numerical calculations 



of Reineltpl 



The intention of the present work, is to offer yet another approach to the fingering 
phenomenon. We will derive a nonlinear differential equation for the finger in the Hele- 



Shaw cell. The equation resembles the Korteweg-de Vries (KdV) soliton equation p 
It will be quite rewarding to find several features emerging from this method, like the 
ability to reproduce stationary finger profiles accurately, and the possibility to deal 
with transient phenomena such as the development and competition between fingers 
and the treatment of tip splitting. In the present work we will focus on the justification 
of the nonlinear soliton-like equation, and its stationary solutions. 

2 Basic equations 

The basic equations governing the fluid flow are described in many works 0. We 
here reproduce them for the sake of completeness, and add our view on the problem. 
The Hele-Shaw creeping flow equation that neglects inertia effects and integrates over 
the transverse dimension of the cell is 01 



-Vp 



12 /i 

where v is the viscous fluid velocity in the plane of the cell, b is the thickness of 
the cell, the viscosity, and p, the pressure in the fluid. Even at this level, there are 
hidden assumptions that are not clearly valid in all situations. The integration over 
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the third dimension indeed helps in reducing the problem to a more manageable set 



of equations. However, the true and real problem is inherently three-dimensional p2 
The existence of films of fluid left behind the finger, puts this aspect in evidence. The 
finger thickness is not constant comparing tail and nose. It is thinner at the nose and 
thicker at the tail. We are therefore already reducing dramatically the complexity, 
and richness, of the phenomenon. The price is paid at the time of comparing the 
measured finger widths as a function of the parameters of the flow to the theoretical 



values. As demonstrated by Reinelt|TB], a careful treatment of films of fluid by means 
of a three (or four) domain splitting of the Hele-Shaw cell, changes dramatically the 
agreement between theory and experiment. 

Originally eq.(|l]) was supposed to imitate Darcy's law for the flow in porous 
mediap|. The phenomenon is nevertheless of a general nature. It belongs to a large 
class of curved front propagation phenomena spanning such diverse topics as dendritic 
solidification and flame combustion. 

The displacing fluid also obeys a Darcian law as above, however, and due to 
the lower viscosity, it is usually ignored and the pressure is assumed constant inside 
the finger. This assumption is also not perfectly justified, but, seems an economical 
working hypothesis, especially when the displacing fluid is a gas. 

The above equation has to be augmented by boundary conditions at the interface. 
Here too the court is open to a variety of conditions, mainly due to the ignorance 
of the third dimension. Basically, the main boundary condition is the one of the 
pressure difference between both sides of the interface. Without corrections of thin 



filmsfTl, 



15|1 the condition reads 



2 T T 

Ap= — cos-f + — (2) 

where, T is the surface tension parameter of the Young-Laplace formula, R, the 
radius of curvature in the plane of the cell and 7, the contact angle between the finger 
and the cell wall in the transverse b direction. Many corrections were proposed to 
the above equation. Park and Homsy[jl5| and Reinelt [jl3| , developed a systematic 



expansion of the effective surface tension parameter as a function of the capillary 
number Ca=^^, where V is the asymptotic velocity of the fluid far ahead of the 
finger that is conventionally related by flux conservation (and without considering 
Poisseuille flow, or else no-slip condition) to the velocity of the flnger U hj V = UX. 
Here A is half the width of the flnger tail after it has fully developed and other 
competing flngers dissolved. 

The term in (||) depending on b is usually dropped, because the flow depends on 
the gradient of the pressure and a constant term does not matter. However, this is 
not valid either. Two effects spoil the constancy, the variation of contact angle 7 
and, the change of flnger thickness as one proceeds from nose to tail, or viceversa. 
It is easy to see that the latter effect is very signiflcant due to the smallness of b as 
compared to the width of the cell. A tiny change in the thickness of the flnger carries 
over to a large modiflcation to the pressure. However, the two-dimensional approach 
can not deal with this matter consistently. 

The contribution of the thinning out of the flnger when coming from the tail to 



the nose adds up to the effect of the surface tension in the plane of the celL One 
can then hope that an effective surface tension might somehow lump up the effect. 
However, this expectation seems to be too optimistic. 



Pitts [|10| found the need to introduce a phenomenological factor to to account 
for the difference in pressure inside the fluid between tail and nose that reflects the 
thinning out of the finger. We will not delve here into the issue, but we do recognize 
that without a reliable treatment of this matter, only the trend of the curve of finger 
width versus capillary number is of interest. In particular we will follow the same 
path that ignores the first term in eq.(0) and use instead an effective (unknown) 
surface tension. For large surface tension, i.e. fingers that cover almost the whole cell 
the treatment will be correct. The narrower the finger as compared to the channel 



width, the stronger the effect of films of fluid left behind [jT3[ Therefore, without a 
three-dimensional treatment we have to be satisfied with results for wide fingers and 
small capillary numbers. Some comments on a possible parametrization that might 
account for the loss of information from the third dimension will be provided later. 
Equation reduces then to 

Ap = I (3) 
with T an effective surface tension parameter. 

To complete the set of equations we assume the fluids to be immiscible. This means 
that there is no interpenetration. The displaced fluid does not cross the boundary. 
Its velocity is tangential to the boundary. The boundary is then along a constant 



stream-line. For a frame of reference with the finger at rest it becomes 



a{v) = tan-^{r]') (4) 

where a is the angle between the fluid velocity and the axis of propagation of the 
finger (x axis) and r]' is the derivative to the finger (y axis) with respect to the x 
axis. The above condition is equivalent to the much simpler one of \E' = constant 
on the finger, where \l/ is the stream function for the viscous fluid p|, , as well as 
to the condition of continuity of the velocity in the normal direction to the interface 



for a moving finger. Note also that the treatment of Reinelt|T^ being structured 
in three different regions, with the front inside an intermediate zone, does allow for 
the crossing of fluid through the front. As it is our aim to provide an alternative 
differential approach, we will not deal with such complications at this stage. 

To the above equations one has to add the continuity equation for an incompress- 
ible fluid of constant density, V ■ v = 0. 

The transformation of the Hele-Shaw cell problem to a two-dimensional situation 
permits the use of analytic functions techniques. Various authors took advantage 
of this feature either to expand the complex flow function in terms of the complex 
coordinate z = x + i y or viceversa0] , or use analyticity conditions to approximate 
the solution [Tl[] The potential $ is defined by v = V$ and the stream function 
by V = Vx\E'. A single component of this vector is needed for flows averaged over 
the cell thickness, "^b, where b denotes the direction transverse to the cell. Moreover 
and in accordance with the Hele-Shaw approach, this function has to be independent 



of the transverse coordinate. The Cauchy-Riemann equations for the potential and 
stream function, imply that the complex potential a; = $ + z "^b, is analytic in z. 

Instead of working with the full complex potential and derive restrictions based on 
analyticity, we here treat the problem by means of the stream function. The stream 
function is a vector, whose curl determines the fluid flow. The incompressibility of 
the flow is then guaranteed. Although the velocity potential $ is a harmonic function 
of the coordinates x, y, as can be readily obtained by applying the incompressibility 
condition to Darcy's law of eq.(|I|), the stream function does not obey the equation 
V^\E' = automatically. The stream function approach and the potential one differ. 
The use of the velocity potential v = V$ assumes a vanishing vorticity, as the curl 
of a gradient vanishes. It is our claim, that this assumption is wrong. The stream 
function approach is superior, as it does not suppose such a neglect of vorticity. The 
velocity potential $ misses a crucial part of the problem. Therefore, it should not be 
considered as a trustworthy tool to find the finger equation of motion. In this respect 
the use of the full potential = $ + i \E'b, is also not valid. Hence, the analyticity 
conditions derived from it implying a harmonic equation for \Ef also, do not reflect the 
full reality of the problem. 

For example, if we consider ^E^^, that depends on x and y only, we have V ■ = 0, 
because \E';, is the z component of the stream vector and it does not depend on z. 
\E'f, satisfies the incompressibility condition, by definition of the velocities. Using 
VxVx^ = VV ■ ^ — V^^, it is clear that, \E'fe will obey the harmonic equation 
if the vorticity, uj = -VxVx^, vanishes. The vorticity may be ignored, except at 
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the interface ||T2| . However, it is there that we need the most the value of the stream 



function. Therefore we will not use the harmonic equation in the present work. 

The boundary conditions alone will suffice for the determination of the stream 
function on the front, and, in turn, allow us to find the evolution equations for the 
finger. However, we will not be able to find the stream function explicitly everywhere 
in the fiuid. 

3 Stream function approach 

The stream function will be determined by resorting to a perturbation expansion. 
The expanded stream function will be made to obey the equations of the fiow as well 
as the boundary conditions. The outcome will be a stream function that is known 
solely at the boundary. This method will enable us to deduce the nonlinear differential 
equations. 

X y 

We first renormalize the axes to x — ^ x = — — , y y = — — . As a consequence 

w/2 w/2 

of this rescaling, the finger transverse dimension becomes less than one, and serves 
as an expansion parameter. 

Two boundary conditions constrain the stream field. The first one is the impen- 
etrability of the cell wall and, the second the immiscibility of the fiuids. The latter 
amounts to considering the interface between the fiuids as a free surface. Dropping 
the suffix b and the tilde, we find 
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9^ 

dx 



, at ?/ = ±1 



(5) 



and 



dt dy dx dx' ^ ^ 

where r]{x,t) denotes the interface curve. Eq.(^ may be easily derived by writing 

drj drj 

the differential of ri as dri = — dt + — dx and, using the definitions 

dt dx 



dx d^ 
dt dy 
dy d^ 
= A = 

drj 

For the stationary finger, in its rest frame we have — = 0. 

dt 

With this substitution, eq.(|0) becomes 



d'^dr] _ _ ^ 

dy dx dx ' 

d^ 9^ r-h 

The differential of \E' d'if = — — dy + —— dx, together with eq. (M) determine the 

dy dx 

unique solution for the stream function in the stationary case to be 

\l'o = constant (9) 
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Without loss of generality, this constant may be taken to be equal to zero, because 
the stream function enters the calculations only through its derivatives. 

For a finger in motion along the x axis, even with a time dependent velocity, but 
still in the creeping flow approximation that neglects kinetic terms, we can write the 
solution to be formally 

^{x,y,t) = ^o{x,y)- J dx ^ (10) 

Eq.(0) solves eq.(|^), with \E'o given by the appropriate solution of eq.(^. 
At the front, and only there, we can rewrite the equation as 

Where f]{x,t) is understood as a function of x and t only. In the general case, 
is no longer a constant on the curved front. 

Consider the stream function in the finger rest frame. The finger is traveling in 
the positive x direction at constant speed U. Moreover, at long distances ahead of 
the finger the fluid is assumed to flow with a constant velocity V. Usually, the no-slip 
condition cannot be imposed in the Hele-Shaw cell and consequently the velocity V 
is considered to be the same all over the cell including at the lateral boundary. This 
could be a good assumption provided the viscosity of the fluid is small, or else the 
lateral edge is not considered as the boundary. However, for large capillary numbers, 
or large viscosities, this procedure seems unreasonable. Boundary layers next to 
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surfaces are patent in common phenomena even for moderately viscous fluids. We 
opt for a more conservative approach and demand the stream line next to the lateral 
edge to carry a velocity V that is not the asymptotic velocity of the fluid ahead 
of the finger. It is a parameter that depends on the flow properties. For capillary 
number tending to infinity it has to be equal to zero. Large capillary numbers may 
be implemented by using very viscous fluids. For such fluids, the no-slip condition 
is a must. Therefore in the limit of infinite capillary number V=0. In the opposite 
limit, i.e. no viscosity at all (finger width going to the full cell width) it is equal to 
the finger velocity. Hence we can assert Q <V <U . Reinelt's equations [0] carry an 
extra term at the lateral edge that allows the velocity V to be variable. 

To our knowledge there have not been any experimental investigations of the fluid 
flow as a function of distance from the cell edge. In light of the previous discussion, 
this measurement could be of relevance. 

Expanding the only component of the stream function, with a separable ansatz, 
antisymmetric in y, so that the finger is symmetric, and in the finger rest frame, we 
find (dropping from now on the tilde on y, i.e. we work with rescaled coordinates) 

^(x, y,t)=y {V{x, t) - U{x, t) + A{x, t) + B{x, t) + C{x, t) + ...) (12) 

with A, B, C unknown functions. We will find below that the expansion will 
effectively be in terms of both y and 1 — y"^. It is then expected to be a reasonable 
approximation in all the cell gap, except perhaps around ?/ ~ 0.5 
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Eq.(^), as well as the constancy of the fluid velocities, for the stationary finger in 
its rest frame imply 



9^ 

Vy = --^ = ,at y = ±1 

V, = —= V-U ,at y = ±1 (13) 
dy 



Inserting eq.([T3|) in eq.(|T^), we are able to fix the functions B and C in terms of 
A. Explicitly 

= A + 3 B +5C 

= A' + B' + C (14) 
primes denoting derivatives with respect to x. 

Recalling that constants are irrelevant for the stream function, the unique solution 
to this order becomes 



B = -2 A 
C = A 



(15) 



The stream function of eq.(0) with the conditions of eqs.(|l^,0) reads 



^(x, y) = {V-U)y + Aix) yil- y'Y + .. 



(16) 



Eq.(|T^) demonstrates that the expansion is a dual one in y and, 1 — y^, as stated 
above. In eq.(^), we have neglected higher order terms and also ignored possible 
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contributions to the stream function from the dependence of the velocities V, U on 
distance. This is clearly not a good approximation for large capillary number as well 
as for the nonstationary problem. 

If we expand the stream function to the next order 

^(x, y,t) = y {V{x, t) - U{x, t) + A{x, t) + B{x, t) + t) + D{x, t) / + ...)(17) 

we need an extra condition to determine the function D. One possible constraint 
to fix this function, may arise from a model of the vortex structure on the curved 
front. The simplest one - although not the most realistic one-, consists in assuming 
a constant vorticity on the interface, with the vorticity defined by = -VxVx^. 
This procedure introduces a new parameter, the vorticity value u. 

In order to stay as low as possible in the number of free parameters, we prefer to 
remain within the expansion to 0{y^) with y < 1. To this order the two conditions of 
eq.d^J^), are sufficient. The constraints we take advantage of, depend on the stream 
function and its first derivative. Higher order expansions as in eq.(p!7D and further 
constraints, as the vorticity, amount to an expansion in second and higher order 
derivatives of the stream function. Due to the fact that, the even polynomials in 
y span a complete set in the rescaled variable < y < 1 for the stationary finger 
case, this expansion presumably converges. The sixth order term is expected to be 
smaller than the fourth order one by factors of the form A^, with A < 1. Although 
the exact dependence in A^ is not known at this stage, it appears a sound working 
hypothesis to truncate the expansion to the sixth order and work with the expression 
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of eq.([l6|). The quality of the finger shapes obtained bellow gives us more confidence 
in this hypothesis. 



4 Stationary and time-dependent nonlinear differ- 
ential equations for the finger 

On the curved front, in the finger rest frame, with the figer nose located at the center 
of coordinates, the condition = constant of eq.(|^) becomes 

^ = 0, aty = r] (18) 
Inserting eq.(|T8D in eq.(|l6D, A{x) is determined at the front to be 



U -V 

^(^) = T, iT2 ' y = Vix) (19) 

(1 -r) 

Using eqs.(|7| J16|JT9|) , the velocity of the fluid at the static finger, the only place 
where we can determine A in closed form, becomes 



Vx,stauc = {V-U)+ A{x) (1 - r/2) (1 - 5 r/^) 
dA{x) 2\2 

Vy, static = r] [l-T] ) 

U -V 

A{x) = -^-^ (20) 

(1 -r) 
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with 7], a function of x. The fluid still flows around the finger in its rest frame, 
hence the name static. 

The rescaled curvature in eq.(0) is 



I v" 

- = r 21 

R (l+r7'2)l 

Where primes denote derivatives with respect to x. 

The sign of the curvature is the appropriate one. This can be seen by using eq.(|l|). 
The left hand side of the equation is a the positive velocity of the finger, chosen here 
to move from left to right, therefore the curvature has to decrease along x. In the 
upper part of the finger the curvature is negative and becomes more so as we proceed 
along X. Therefore the right hand side of eq.([l|) is positive with the positive sign in 



eq-dH). 

Eq. (^) is the standard expression for the curvature, and, it is easily derivable from 
the definition of the arclength and the corresponding angle for differential increments. 
Equations (||,^^1]) after rescaling imply 



fb^dl 



= - 



T]" 



I2\ 



(22) 



3 w'^n dx 

The stationary finger equation is obtained now by evaluating the derivative in 
eq. (|2^) , with 



4 rf 

Vx.static = (y ~ U) 

1 — 



(23) 
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from eq.(^); and transforming back to the rest frame of the cell, with the finger 
in motion, = v^^static + U. 

After some straightforward algebraic manipulations we find 



4 e + 1 — 5 r]^ 

1-7^2 



W{v) = 4^-^ (24) 



where the suffix indicates differentiation with respect to x. 

In terms of the arclength measured from the tail of the finger, the equation be- 



comes 



I — rjg 4 B 

(25) 

V 

In eqs.(0,^), e = — . As was found in previous works 0, || the only parameter 

1 . . 1 12/i U TT 1 Ol 

entermg the equation is — = with w being the width of the Hele-bhaw 

B T h"^ 

cell and all the variables are scaled to units of half this length. The existence of a 
single parameter is clearly insufficient as demonstrated by Reinelt [|T^ . Customarily e 
of eqs. (|2^j25|) is equated to A, half the finger width. We here take it as a parameter as 



explained in section (j^). Asymptotically far back at the tail we must have W{ri) = 0, 
5 

4 A2 



or, e = — — . The parameter e, is then determined by the solutions to the 



equations and is not a phenomenological parameter. 
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The potential W guides the propagation of the finger. The equation is of third or- 



der in the spatial derivatives as the Korteweg-deVries equation [|14[. Third order differ- 
ential equations for interface propagation are well known in the literature. Some noto- 
rious examples are: Landau and Levich[0], Bretherton|[T^, and Park and HomsyJlS 



The third order equation (|^) can be transformed to a more manageable second 
order one in terms of the angle tangent to the curve, with rj obtained by integration. 
For ds starting at the tail, where 9 ^ n, ds cos{6) = — dx. 

Equation (EE) now reads 



= ^-co.w'^'f"' 



(9s2 ^ ' AB 
X- j ds sin{e) (26) 



With W{ri) defined in eq. 



The time dependent equation is obtained by using eq.(|Tl|), and recalling the def- 
inition of Vx in terms of the stream function of eq.(|^. For the velocity in the x 
direction, we now have 



Vx = -— + U + Vx,static (27) 



w 

Rescaling as before the length coordinates t], and x by — , and the time by 
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t 



t u 



w/2' 
equation 



and inserting eq.(^) in eq.(^), we obtain the time dependent differential 



3^ + {I + Vl 



4 B 



Wir]) - ^ 



In terms of the arclength eq.(p8D becomes 



r]sss + -, 3 + (l-r7j 



4 B 



W{r]) 



\ 2^* 



(2J 



(29) 



While for the angle tangent to the front, the time-dependent equation ( ]28|) reads 



5?-^"^^^) 45 
\ — ds sin{9) 



wi,) + cosie) -/- 



(30) 



To gain confidence in the equations, we proceed to show that the time dependent 
equation (^) gives the correct dispersion relation for the perturbation of a planar 
front. We will find that the instability parameter in time, here called a, will obey a 
dispersion relation of the form a = \a\ (1 — Aa^), with a related to the wavenumber of 
the perturbation and A a parameter proportional to the surface tension. The velocity, 
or pressure gradient, destabilizes and the surface tension stabilizes. We will connect 



to Chuoke et al.[l£l in the derivation. 



We will consider a perturbation protruding from a moving front with velocity U. 
Eq.(p8D is firstly rewritten around y ^ 
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This equation is equivalent to 



dx 



X = 4 5—^ +x - Idx^ (32) 



2 



Only around y ~ 0, the equation can be integrated once analytically. 



Now we need an ansatz for the solution. Following Chuoke et al.|T9[ we use the 
expression 



x{y,t) = Pia)e'f' 



a y + a t 



(33) 



(5 is the amplitude of the perturbation, that is taken as dependent on a due to the 



nonlinearity of the equation, as Chuoke et al.[jT9l assumed. (See below their equation 
[4]). 



The integral in eq.(|3^ may be rewritten as 



/ = / dx^ 

Vx 



Xy Xt dy at y = r] 



(34) 



Evaluating the integral (B3) with eq. (p3|) we find 
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2 

a P'^a y 



(35) 



The curvature may be calculated now by using 



X. 



yy 



-X. 



yy 



-c? X 



(36) 



Using eq. (|36D , the integral of eq.(p5D, and, x of eq. (|33|) , at t = to 0{y ), as we 



are assuming a perturbative expansion, x oi eq.(32), becomes 



Xo 



-4 B a'^p + p) {1+ ay) - a p^a y 



(37) 



The constant piece of xo in eq. (|37|) , —4 B a^P + /?, is irrelevant. This is the 
reason we do not need to specify the lower bound of the integral in eq.(0). We took 
it as y = 0. Any other choice will merely change the unimportant constant in x- 



To satisfy eq.(p2|), to the lowest order in j/, eq. (|37|) has to obey the algebraic 
condition 



AB a^P + pa- ap^a = 



(3J 
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The third term in eq. (|38|) is quadratic in /3. Independence from the initial pertur- 



bation amphtude, requires /5 to be a function of a, as assumed in eq. (|33D . We can 
now determine this dependence to be Pa = k, with k, a constant. /3 is a positive 
number, therefore, k has to be positive for positive a, or negative for negative a. 
Equivalently, the proportionahty above can be written as 



P = ri- (39) 



\a\ 



Inserting eq. (|39D into eq. (|38D we find 



\k\ 

-AB a^l3 + /5a - = (40) 



\a\ 



Finally, redefining lal ^ t-t, A = 4 5 /c^ we find the expected dispersion relation 

\k\ 



a=\a\{l-\ a^) (41) 



The time- dependent equation was developed on the basis of an antisymmetric 
stream function, eq. (p!6D , the solutions, may, nevertheless, lack any predetermined 
symmetry. This phenomenon may be related to the physical concept of hidden sym- 



metries and spontaneous symmetry breaking|20|. Although the equations are sym- 
metric under a certain transformation, the solutions can break the symmetry. Bound- 
ary conditions and initial conditions are key factors in determining the symmetries 
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of the solutions. Such is the case of an initial asymmetric perturbation to the front. 
On the experimental side, a symmetric finger develops only asymptotically at long 
times. The transient behavior is patently asymmetric. 



5 Stationary finger solutions 

Let us analyze the predicted finger half-widths A, as a function of the parameter B. 

As discussed in section (|), we expect the parameter e to be equal to one in the limit 

of -B — > oo, that is possible to implement with a displaced fluid of zero viscosity, 

but a density that still assures the instability and creation of a finger. In this limit 

the potential determines the finger half-width ( location at which the derivatives of 

rj vanish) to be A = 1, as expected. On the other end of -B ^ 0, or infinite viscosity 

limit, we argued that e = in the spirit of the no-slip condition for viscous fiuids. In 

this limit A = ^ ~ 0.447. Looking at the experimental work of Tabeling et al.0, 

especially figure 8, it is seen that for — ^ oo the finger width for various aspect 

B 

ratios tends to approximately A = 0.45. 

Therefore we can be quite confident on the potential W. On the other hand, the 
dynamics as described by the third order differential equation linear in the parameter 
B is certainly only approximately correct. In the same manner as the KdV equation is 
only an approximate equation derived by expanding to a certain order in the deriva- 
tives. The need of higher order terms in the curvature is apparent from the work of 



Brower et al.p2[. As stressed in the introduction, a full three dimensional account of 
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the fluid flow is necessary in order to take care of the thinning out of the finger in the 

transverse direction. Alternatively we could view the equation as essentially correct, 

with the dependence on B nonlinear and unknown. Consequently, we can trust the 

equation in the limit of small — . 

B 

The numerical solution was started from positive y = X and followed around the 
finger to y = —A. We employed a fourth order Runge-Kutta algorithm. We opted 
for the mixed integrodifferential approach of eq. (p6D , that proved easier to handle 
numerically. The equation can be solved without this change of variables also with 
the full third order equation (P^. 

For fixed B we varied manually the value of the asymptotic half-width (thereby 
fixing e), the slope 6' and the value of 6 less than 6 = 180°. The equation is integrated 
both backwards and forwards and the parameters are fine tuned until the shape of 
the finger looks as symmetric as possible x{—X) ~ x{X), and the nose tip shows up 
at exactly 90'^. We inforced the symmetry condition to an accuracy of around 1%. 
Moreover, only solutions where the angle decreases smoothly from the tail forwards 
were accepted. We rejected meandering (even slightly so) solutions that could appear 
as a many prongued finger. 

This method of tedious search proved to be reliable. Other algorithms, such as a 
shooting method, etc., ended up in unreasonable, asymmetric solutions. 

The literature does not abound in theoretically calculated profiles of fingers. We 
could hardly find any such solutions besides those of McLean and SaffmanR], the phe- 



nomenological solutions of Pitts|TO|, and, the original ones of Saffman and Taylor [| 
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Figure 1: Rescaled Gnger proGle for X = 0.93 at = 12. Numerical solution, solid 
line, analytical solution of Pitts [10], dashed line. 

The latter obtained without the inclusion of surface tension effects. It appears then, 
that it is not an easy task to predict properly finger shapes. 

We proceed to show 3 profiles found for finger half-widths A = 0.93, 0.82, 0.66. 
They are compared with the analytical solution of Pitts [|10| that fits the data ex- 
tremely well for finger half-widths smaller than around A ~ 0.8 and lies slightly 
below the data points for higher values of A. The figures are depicted for rescaled 
coordinates in terms of A as has became customary 0, |T0| . 
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Figure 2: Rescaled Anger profile for A = 0.82 at ^ — 20. Numerical solution, solid 
line, analytical solution of Pitts [10], dashed line. 
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From figures 1-3 it is clear that the solutions match the phenomenological solutions 
of Pitts cos(— — )e~2~ = 1, and are even superior to it for large values of A. The 

2 A 

dependence of the asymptotic finger half-width A on B, is shown in figure 4. The 
present results improve upon the calculations of McLean and Saffman]^ in the range 
of B for which we were able to find reliable solutions, but still lie below the measured 
values that were obtained from averaging the results of Tabeling et al. . for various 
aspect ratios. Note however, that the results do not include finite gap corrections that 
lift the curve upwards, whereas the results of McLean and Saffman]^. shown in the 
figure do adjust for such corrections. Without these corrections the latter results fall 
much below the depicted curve. As stressed above, finite film effects are extremely 
important, as evidenced from the calculations made Reinelt | T3[ . 



We found numerically stationary finger profiles for B < 110. Beyond this value, 
the large contributions introduced by the potential in eq.(p6|), render the integration 



from the tail unstable. This is not unexpected in light of the arguments provided 

in section (0). The curve of A versus — for large values of this parameter, can be 

B 

fitted by a function of the form B e'' ^ , with k a parameter. In future investigations 
higher order terms in the potential will be added, in the search of such a functional 
dependence in a self-consistent stream function approach. 

Several attempts were made in the present study, to improve the algorithm by 
matching solutions run from the tail and from the nose. However, starting from the 
nose is very problematic, because of the singularity in the slope, that no matter how, 
creeps up in any scheme, even after transforming to better behaved variables. 
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Figure 3: Rescaled Gnger profile for X — 0.66 at ^ — 48. Numerical solution, solid 
line, analytical solution of Pitts [10], dashed line. 
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Figure 4: Asymptotic half-width of the finger, X as a function of -g. Present work, 
full line, theoretical results of McLean and Saffman [7] corrected for finite film ef- 
fects, dashed line, and, experimental results of Tabeling et al. [9], dash-dot line with 
squares. 
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The equation is extremely sensitive to the slope and initial angle that must be 
close to 180°. Fortunately, it was found that for a a fixed value of B, there cor- 
responds a single value of A; a unique finger solution. The selection of the finger 
half-width is univocous, with all the numerical methods employed. However, we can- 
not rule out completely the existence of other unstable branches not accessed by the 
numerical method, as found by Vanden-Broek|p for the McLean and Saffman[0 set 
of integrodifferential equations. 

In the course of the numerical investigation we found solutions that correspond 
to various multiple finger-like structures. These structures originate from specific 
boundary conditions at y = \. The existence of such complicated patterns is quite 
encouraging, especially if we want to proceed further to the production, development, 
competition and splitting between fingers with the time dependent equation displayed 
above0. 



6 Conclusions 

The contribution of the present work consists in the development of nonlinear differ- 
ential equations for the Hele-Shaw Saffman- Taylor fingering problem. The approach 
is based upon the treatment of the stream function as a separable potential without 
any free phenomenological, adjustable parameters. The numerical stationary profiles, 

match extremely well the experimental ones, that are, in turn, quite accurately re- 
^Beautiful martingales similar to the twisting and whirling of plant branches were also found. 

See remarks on symmetry below eq.(p8|) 
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produced by the phenomenological solution of PittsfT^. The predicted dependence 



of the asymptotic finger half-width A on the physical parameters of the problem, is 
also satisfactory. 



For the stationary case, the finger nonlinear equations are eqs.(]2^), and (^) the 
latter expressed in terms of the arclength, while eq. (p6D is the corresponding integrod- 
ifferential equation, of second order, for the angle as a function of arclength. These 
are nonlinear and nonlocal equations. 

The time-dependent equation is of the third order in space, and first order in time. 
Eqs.(P^,P^), corresponds to t] in terms of x and the arclength respectively, whereas 
eq. (|30|) is its integrodifferential equation for the angle as a function of the arclength 
also. 

The equations found here, resemble the Korteweg-deVries (KdV) equation, that 



has been long recognized as a solitary wave equationp3[. This equation was obtained 
for the nonlinear propagation of shallow water waves in a channel , and it is of the 
same order in space and time, as our equations. 
The KdV equation reads 



fh^ 2Th\ 3 Co , rit^ r. 

'^0 T ]Vxxx + ttt'^Vx + Vx{co H ) = (42) 

V 6 p g J 2 h 7]^ 

where h, is the depth of the channel in which the soliton propagates. Cq = y/g~h 
is the speed of linear propagation of waves in the channel, g is the acceleration of 
gravity and rj denotes the soliton profile propagating in the channel along the x 
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axis with velocity U = , in a fluid with density p and surface tension T. 

Vx 

This is a third order equation nonlinear in t]. The specific nonlinear term, and, 



the potential found in eq.(^J), differ from the quadratic nonlinear term of the KdV 



equation (^2|). The nonlinear terms in eq.(^), and eq.(^) arise from the boundary 
conditions at the interface between the front and the displaced fluid as well as the 
surface tension pressure drop. Eq.(p^ however, carries also the information of the 



lateral flow constraint that is absent in eq. (^2]) . Moreover, eq.(^) is based upon the 
creeping flow assumption of eq. (|^), while eq.(^2D does not assume such a restricted 
version of the Navier-Stokes equation. Therefore the nonlinear terms differ. 

Nevertheless, the equations are similar. Both are of third order in the shape 
variable, both are nonlinear, and are linear in time In the derivation of eq.(P^ we 



adhered to the separable assumption expansion used by Korteweg and deVries||23| and 
their techniques. It is then, perhaps, not surprising, that we found a similar equation 
- despite the differences in context-, as both equations describe the propagation of 
stationary and stable shape fronts in a fluid medium. 

Links between soliton equations and the Hele-Shaw finger are already hinted in 
the work of Kadanoff[|T^, in which varieties of the so-called Harry-Dym equations 
were found to be in close relationship to finger development. The above equations 
belong to a broad class of Eux-like partial differential equations, that are widely used 
in the literature. In the context of Darcian flow in the Hele-Shaw cell, Goldstein et 
al.[pT| study instabilities and singularities by means of an equation of this type. The 
main difference between the equations we derived here, and the one of Goldstein et 
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aLpT], is, again, the appearance of a potential term W{ri) in eqs.(p8t 

The present method does not use any phenomenologically observed characteristic. 

The model results follow solely from a treatment in terms of the stream function, and 

the relaxation of the condition of A = that amounts to an unphysical restriction on 

the Hele-Shaw cell walls especially for large viscosities, or small capillary numbers. A 
5 

4 



flat tail requires e = — — . e, is then determined by A, that is in turn determined 



by the existence of solutions. We assumed a constant, capillary number dependent e. 
Further progress requires a model for the dependence of e on the capillary number. 

The predicted dependence of finger half-width A, on capillary number through 
the parameter B of figure 4, achieves more than the integrodifferential approach of 
McLean and SafFman^ even without corrections of finite films left behind the prop- 
agating finger. The partial success in dealing with the problem, suggests that front 
propagation phenomena might have a corresponding soliton-like guiding equation, 
at least for restricted geometries. 

For the Hele-Shaw cell, the cell provides the potential for the propagation of the 
finger, while the pressure difference determined by the surface tension provides the 
dynamics. In dendritic growth, the temperature gradient provides the dynamical 
guidance and the cell and substrates, the potential. 

The obvious shortcoming of the current method, and other existing methods, is 
the neglect of the dynamics along the third dimension. One possible way to improve 
the equations would be to include vorticity. Even if in the bulk there is no vorticity, 
on the front there could be a sizable amount of it without spoiling the bulk properties 
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of the fluid in the hydrodynamic regime, as in superfluid HehumP4|. It seems then, 
that a treatment including vorticity may be more appropriate. 
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